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By fine-tuning generic Cauchy data, critical phenomena have recently been discovered in the 
black hole/no black hole “phase transition” of various gravitating systems. For the spherisymmetric 
scalar field system, we find the “critical” spacetime separating the two phases by demanding discrete 
scale-invariance, analyticity, and an additional reflection-type symmetry. The resulting nonlinear 
hyperbolic boundary value problem, with the rescaling factor A as the eigenvalue, is solved numer¬ 
ically by relaxation. We find A = 3.4439 ± 0.0004. 
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Recently, Choptuik [4] has studied the gravitational collapse of a real scalar field (massless or massive, minimally 
or conformally coupled) in spherical symmetry, using an adaptive mesh refinement numerical technique which allows 
him to study details on very small spacetime scales. To describe his results concisely, we invoke coordinates {p,p} on 
the phase space of the spherisymmetric gravitating scalar field, where p is any smooth coordinate such that p = 0 is 
the hypersurface which divides black-hole from no-black-hole spacetimes, while p denotes the remaining coordinates. 
Choptuik’s results strongly indicate the following conjectures: 

(1) For any choice of coordinate system {p,p}, the mass of sufficiently small black holes is given by M = f(p)p' y 
(“scaling”), where 7 ~ 0.37 is a universal exponent. 

(2) There is a “critical solution” {p = 0 ,p = p*(t)}, which acts as an intermediate attractor in a thin sheet 
surrounding the p = 0 hypersurface on both sides (“universality”). 

(3) This solution shows a discrete homotheticity (“echoing”), or scale invariance, to be defined more precisely below, 
with a logarithmic rescaling factor A ~ 3.44. 

More recent research indicates that these properties hold for other self-gravitating systems. Universality, echoing 
with (A ~ 0.6) and scaling (with 7 ~ 0.37) were found in collapse of axisymmetric gravitational waves 0 . Universality, 
continuous self-similarity and scaling (with 7 ~ 0.36) were found in perfect fluid collapse with p = p/3 [g). The exactly 
self-similar solution was calculated as an eigenvalue problem ||, and the critical exponent calculated to high precision 
(7 = 0.3558019) by perturbing it Q. Critical exponents for other values of k in p = kp were calculated in this 
manner j^] and confirmed in collapse calculations ji|, with values strongly dependent on k. Self-similar solutions were 
also calculated for a complex scalar field [[tJ and an axion-dilaton combination Q. A value (7 = 0.387106) of the 
critical exponent for the complex scalar field was derived by perturbing the self-similar solution |J , but the latter is 
apparently not an attractor ©• 

Universality, scale invariance and critical exponents indicate an exciting new connection between renormalisation 
group theory and classical general relativity. 7 appears to vary from one physical system to another, but its values 
for vacuum or trace-free matter are remarkably similar in the examples found so far. 

In this Letter, we impose echoing, and an additional reflection-type symmetry, in our ansatz, together with analyt¬ 
icity, and solve the resulting nonlinear hyperbolic eigenproblem, instead of evolving and fine-tuning Cauchy data. In 
the language of renormalisation group theory, we find a fixed point of gravitational collapse under a rescaling of space 
and time by solving the renormalisation group equations. In a future paper we intend to calculate 7 by perturbing 
around the fixed point, along the lines of y|||J. 

The Einstein equations we consider here are 


G a b — 87 tG I (j),a/),b 0 {Jab/),erf*’ 


Cl) 

in spherical symmetry. The matter equation 4> c ' c = 0 follows as a Bianclii identity. Following Choptuik, we define 
the metric as 

ds 2 = —a(r, t) 2 dt 2 + a(r , t) 2 dr 2 + r 2 ( dO 2 + sin 2 0 dp 2 ), (2) 

where the remaining gauge freedom is fixed by the condition a(r = 0, f) = 1 , and auxiliary matter fields as 
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( 3 ) 


X(r, t) = \FZkG — <fi r , Y(r, t) = V2nG — (f> t . 

a ’ a ’ 

The symmetry of the attractor observed by Choptuik can be expressed in coordinate language as Z(r,t) = 
Z(re A ,te A ), where Z stands for any one of a, a, X and Y, and A ~ 3.44 is a constant. Here the zero of t has 

been adjusted so that (0,0) is the accumulation point of the echos. We introduce auxiliary (nonmetric) coordinates 

where the symmetry appears as a simple periodicity of the Z: 

r=ln(±t), £ = ln (±r/t), Z(£,t) = Z(£,t + A). (4) 

In the following, we use the upper sign (t > 0). The equations for the lower sign are obtained by changing the sign 
of Y in all following equations. 

Geometrically, the symmetry can be described as a discrete homotheticity: when we Lie-drag g ab along the vector 
field d/dr by the distance A, we obtain g ab e 2A , while T ab is mapped to T ab . The vector field along which we Lie-drag 
is not unique, however, because we only ever consider the effect of Lie-dragging the finite distance A. We parameterize 
this arbitrariness by introducing a free periodic function £o ( r ) into the coordinate system such that the vector involved 
in the symmetry is still d/dr. At the same time, for clarity of presentation, we absorb A into the coordinate r. We 
therefore define the coordinates in which we are going to work as 

p = 2nr/A, C = £-£o (t), Z((,ip + 2tt) = Z((,ip). (5) 

Evans and Coleman |j| found the critical spacetime of spherical fluid collapse by imposing continuous homotheticity. 
In our coordinates, a, a, and the fluid variables corresponding to X and Y are then functions of £ alone, the Einstein 
and matter equations are reduced to a system of nonlinear ODEs, and the solution is uniquely specified by regularity 
conditions | 0 . Here we use a similar approach in order to find Choptuik’s critical spacetime of scalar field collapse. 
(After this work was begun, the same approach was used to calculate (continuously) self-similar solutions for the 
complex scalar field |Q and an axion-dilaton combination ||.) 

We solve the field equations for the ^-derivatives of the fields Z as functions of these fields and their (^-derivatives. 
It is convenient to use the new field g = e^^a/a instead of a , and X± = X ± Y instead of X and Y. The resulting 
equations are 


o,c = [(! - « 2 ) + a 2 (X\ + X 2 _)\ , 

g, c = g( i-« 2 ), 


*+,c 

X-,c 


B+ 

1 + D' 
B_ 

1 -D 1 


where we have introduced the abbreviations 


z= 1 + 


27t d £o 

A dip 


B±= -(1 — a 2 )X± 


, D = z 1 e c 5 , 
a 2 XlX± - X T ± z ^DX ±!ip . 


( 6 ) 

( 7 ) 

( 8 ) 
(9) 


( 10 ) 

( 11 ) 


There is also one equation containing only the Z and Z tip , 

z ^ a ~f = \ [(1 - a 2 ) + ^ {*% + X -) + ( X l - X -)) ■ ( 12 ) 

It acts as a constraint, which is conserved by the four “evolution equations” above. 

For small enough ( these equations define a constrained Cauchy problem, with ( playing the role of time, on 
the cylinder obtained by identifying ip with period 27r. At £ = — oo, corresponding to r = 0, we set the boundary 
conditions a = 1 (regularity of the metric) and a = 1 (coordinate condition). Expanding the field equations in powers 
of e^, we find that data obeying these conditions are determined by fo^) and one more free function Yo^): which is 
defined by the expansion 

Y{<p, C) = Y 0 (ip)e ^ e f + O (e 3C ) . (13) 
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As £ increases, the Cauchy problem eventually becomes degenerate, when D = 1. In analogy to the “sonic point” 
of the ODEs describing continuously homothetic spacetimes in , we call this line the “sonic line”. The equation of 
radial null geodesics is 


dc 

dp 


A 

27 TZ 


(-1 ±D-1). 


(14) 


The sonic line is therefore the set of points where a null geodesic touches a surface of constant £. (In general, it would 
be a matter characteristic, but for our choice of matter these are identical with the null geodesics.) The solution can 
be uniquely continued across the sonic line when we impose analyticity. As a technical simplification, we make use of 
the coordinate freedom in £ 0 by moving the sonic line to £ = 0. Then we can enforce analyticity simply by expanding 
in powers of £. We find that regular data near £ = 0 can be expressed in terms of £o {<p) and one more free function 
X + q( p), which is defined as 


X+o (<p) = *(£ = 0 ,p) + Y( £ = 0 , p). (15) 

We have now formulated a hyperbolic boundary value problem on a rectangle with two sides identified (a finite 
cylinder) in 1 +1 dimensions. We have three independent fields, for example X +1 X_ and g. On the other hand there 
are three free functions in the boundary data, minus one degree of freedom corresponding to translations in p, plus 
A as the eigenvalue of the problem. By this count we expect solutions to be locally unique, with a discrete spectrum. 

We cut the number of degrees of freedom in half by imposing the additional symmetry Z(p + n) = ±Z(p), with the 
+ sign holding for a, g and £o, and the — sign for X+ and X-. It is consistent with the field equations and Choptuik’s 
data [jT2|. Moreover, Choptuik observed that the massive scalar field has the same attractor as the massless one 
considered here. The necessary and sufficient condition for this is that <f> remains bounded, because the mass term in 
its stress tensor is then dominated by the gradient squared term as <f> varies on ever smaller spacetime scales. For <f> to 
remain bounded its derivatives X and Y must have vanishing zero frequency Fourier components in p, which in turn 
requires that all their even frequencies vanish, or else these could be combined to give a zero frequency contribution 
in the evolution equations. It follows in turn that a, g and £o must not have odd frequency components. 

As we are dealing with smooth periodic functions, it is useful to decompose all fields into their Fourier components 
with respect to p. Integration and differentiation are done in Fourier components. Algebraic operations are done in 
p space, which makes our algorithm pseudo-spectral. Due to the nonlinearity of the problem, dealiasing turns out to 
be essential for stability. We dealiase convolution sums by using a number of collocation points in p equal to twice 
the number of Fourier components. 

Because the number of variables is large and the problem is nonlinear in an essential way (there are no regular 
solutions to the linearized, no gravity, problem), any algorithm is likely to have only a small region of convergence 
around any solution. We have therefore started with an initial guess sufficiently close to Choptuik’s critical spacetime, 
in order to establish that this solution exists having the echoing symmetry as an an exact symmetry, and that it is 
locally unique, and to calculate it with higher precision than has been possible by fine-tuning Cauchy data. A global 
search for all solutions that may exist is desirable, but not possible with the present algorithm. 

Our manual input into the algorithm is limited to the following guess for lo and A: lo = — 2.3 sin — 0.6sin3<p, and 
A = 3.44. Here and in the following we fix the translation invariance by defining lo to have no cos<p component. The 
numbers were estimated from very near-critical collapse data made available by Choptuik | jl2| , after transformation 
to coordinates (£,r). In a first step, we begin with the very rough guess £o = 0, and shoot from £ = —oo towards 
increasing £. When D(p, £) first gets close to 1 in two points p, and is therefore about to become singular in 
those points, we stop the evolution and calculate a new value of £o that is designed to “flatten” D{p), i.e. to make it 
roughly D(p) ~ 1 for all p at that £. Then we shoot again, thus iteratively improving £o- After convergence, we read 
off A' + at the endpoint of our one-sided shooting, which by now is close to £ = 0, and thus have an initial guess for 
A +0 as well. 

As an intermediate step, we calculate an initial guess for the values of all fields on a grid in £ by shooting from 
both £ = —oo and £ = 0 to a fitting point, typically £ = —1. This involves a Taylor expansion around the regular 
singular point 0 as well as around —oo. Using this expansion and shooting from £ = 0 transfers the bulk of the error 
in our improving solution away from the point £ = 0 to the fitting point £ = — 1 , making it easier for the following 
step to handle. 

In the last step, we go over to a standard relaxation algorithm (|l3|. For the purpose of relaxation the independent 
variables at each grid point in £ are the odd Fourier components of X + and A_ and the even components of g and 
£o- a is not considered as independent, but reconstructed at each step from the other fields by solving equation dis¬ 
solution of this ODE is by iteration of the corresponding integral equation. (The constant component of a has to be 
calculated separately.) Between generic grid points in £ we enforce the discretized £-derivatives 
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A. 


n+1 


- Z n = hF 


-"n+l 


+ z„ 


(16) 


(The ([-derivative of £o is zero by definition). At the boundary £ = —oo we enforce relations between g, X, and Y 
derived from expanding the field equations, and at £ = 0 we enforce D = 1 and B_ = 0. The relaxation part of the 
algorithm is much simpler than the shooting parts, but the latter appear to have a larger region of convergence, thus 
serving as a stepping stone. 

The boundary data of the solution have been tabulated in Table 1. In particular, the echoing period is A = 
(3.4439 ± 0.0004). The error bars have been obtained combining the results of three different convergence tests. 
(1) We compare the results obtained for different numbers M of grid points in (. As expected the convergence is 
quadratic, over a wide range of M, but only up to some maximal value. (2) The convergence of the tabulated numbers 
with increasing number N of Fourier components used in the calculation is rapid (“spectral convergence”) for N > 32. 
(3) £ = — oo is represented by a finite value of £, using a Taylor expansion to one beyond leading order in exp£. As 
expected, this convergence is quadratic in exp £, over some range of C- As long as the difference between runs of 
different precision has the expected functional form, we can use it to estimate the numerical error. The tabulated 
data are from a run with M = 201 equally spaced points in the interval —5 < £ < 0 and N = 64 components (half 
of which vanish) per function, compared with —6 < C < 0, M = 401 and N = 128 respectively. The three sources of 
numerical error are comparable for this choice. 

We have compared the fields a, a, X and Y with Choptuik’s data, after interpolating to the largest rectangular 
grid in r and £ contained in both data sets, with —3.2 < £ < 1.3. We have evaluated the root mean square of the 
absolute difference point by point of the fields a, A' and Y (which are bounded and of order one in the solution) and 
the relative difference in a (which is unbounded above, but bounded below by 1). After adjusting a non-universal 
offset in r between the data sets, this difference is 3.9 1CU 2 for a, and somewhat smaller for the other fields. By 
comparison, the estimated root mean square pointwise error in our data is 1.6 10 -3 in a, and 1.0 10 -4 or less for the 
other fields. We have therefore improved the precision with which the Choptuik spacetime is known by one to two 
orders of magnitude, while A is now known to one part in 10 4 . Future improvements are possible. 


Data files of the solution are available at ittp://www. laeff . esa.es/~gundlach/ 
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TABLE I. Decomposition of the boundary data in sines and cosines. Their period is A = (3.4439 ± 0.0004). 


Component 


Component 

Vo 

A+o 

constant 

(1.5813 ± 0.0007) 

COS f 

0 a 

(-4.3831 ±0.0006) lO^ 1 

cos 2ip 

(6.658 ± 0.006) 10' 2 

sin f 

(-2.364 ± 0.006) 

(-3.2287 ± 0.0008) 10" 1 

sin 2 ip 

(-1.577 ± 0.002) 10" 1 

cos 3 ip 

(-1.46 ±0.05) lO^ 1 

(6.74 ±0.05) 10“ 3 

cos 4 ip 

(-2.014 ± 0.004) 10” 2 

sin 3(p 

(-9.52 ±0.08) lO^ 1 

(1.017 ±0.001) 10" 1 

sin 4y> 

(-3.3 ±0.2) 10“ 4 

cos 5 ip 

(—1.12 ±0.05) lO^ 1 

(2.431 ± 0.003) 10~ 2 

cos 6 ip 

(1.979 ± 0.005) 10“ 3 

sin 5 tp 

(-4.06 ±0.06) lO^ 1 

(-1.807 ±0.008) 10~ 2 

sin 6 ip 

(2.249 ± 0.008) 10“ 3 

cos 7 ip 

(-7.0 ±0.4) 10~ 2 

(-9.85 ±0.02) 10“ 3 

cos 8 ip 

(1.37 ±0.01) 10“ 4 

sin 7 ip 

(-1.73 ±0.04) lO^ 1 

(-2.14 ±0.02) 10“ 3 

sin 8 ip 

(-8.186 ± 0.004) 10“ 4 

cos 9 ip 

(-3.9 ±0.3) 10~ 2 

(1.76 ±0.01) 10“ 3 

cos 10 ip 

(-1.886 ± 0.002) 10“ 4 

sin 

(-7.3 ±0.2) 10~ 2 

(3.116 ±0.005) 10“ 3 

sin 10 (fi 

(3.6 ± 0.1) 10 -5 

COS lllfi 

(-2.1 ±0.2) 10~ 2 

(3.97 ±0.02) 10“ 4 

cos 12 ip 

(3.08 ±0.03) 10“ 5 

sin Hip 

(-3.0 ±0.1) 10~ 2 

(-1.195 ±0.005) 10“ 3 

sin 12ip 

(1.73 ±0.02) 10“ s 

cos 13ip 

(-1.1 ±0.1) 10~ 2 

(-4.05 ±0.01) 10“ 4 

cos 14 <p 

(-4.05 ±0.09) 10" 6 

sin 13y 

(-1.24 ±0.06) 10~ 2 

(1.93 ±0.02) 10“ 4 

sin 14ip 

(-1.371 ± 0.002) 10“ 5 

cos 15ip 

(-5.2 ±0.9) 10~ 3 

(1.53 ±0.01) 10“ 4 

cos 16 ip 

(-3.05 ±0.03) 10~ 6 

sin 15ip 

(-5.0 ±0.2) 10^ 3 

(5.94 ±0.06) 10“ 5 


a by definition, to fix translation degree of freedom 
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